//@+leo-ver=5-thin
//@+node:sonam.20180318222424.1: * @file SysElmnt.pas
//@+at
//  該單元實現了體繫派的五要素計算
// 
//@@c
unit SysElmnt;

//@+others
//@+node:sonam.20180328221954.1: ** Interface
interface

uses SysUtils, FracData, PlanetCycle, SunMoonTab;

Type
//@+node:sonam.20180328222143.1: *3* 類聲明
	EZeroYearException = class(Exception);
	TSysElement = class(TObject)
    public
		constructor Create(bk: CalcBook);

		function AccumMonth() : Integer;    // 積月
		function MonthRemainder() : Integer; // 閏餘
        function IsExcessMonth() : Boolean;   // 是閏月嗎?
		function PlanetBase() : TWeekData; // 曜基數
		function IntFraction() : TMoonData;  // 整零
		function SunBase() : TSunData;      // 太陽基數
		function PlanetMid() : TWeekData; // 中曜
        function SunMid() : TSunData;       // 中日
        function MoonStep() : TWeekData;  // 月步
        function SemiAccuratePlanet() : TWeekData; // 半定曜
        function SunStep() : TSunData;      // 日步
        function AccuratePlanet() : TWeekData; // 定曜
        function AccurateSun() : TSunData;    // 定日
		function MoonStar() : TSunData; // 太隂日月宿
        function SonarMoonStar() : TWeekData; // 太陽日月宿
        function Intersect() : TSunData; // 會合
        function Effect() : Integer;    // 作用
    private
		fBook : CalcBook;

		fYear : Integer;
		fMonth : Integer;
		fDay : Integer;
        fLeapMonth : Boolean;
		fCycle : TPlanetCycle;
		MoonBalance : Boolean;
		SunSecondHalf : Boolean;
		procedure SetYear(y: Integer);
		procedure SetMonth(m : Integer);
		procedure SetDay(d: Integer);
        procedure SetLeapMonth(l: Boolean);
        function PlanetAccumMonth(y: Integer) : Int64;
	public
        property Year : Integer read fYear write SetYear;
        property Month: Integer read fMonth write SetMonth;
		property Day : Integer read fDay write SetDay;
        property LeapMonth : Boolean read fLeapMonth write SetLeapMonth;

        //此處的y雖然有勝生年、西元紀年和自丁卯已過週年，m和d都是藏歷月和日
		procedure SetDate(y, m, d: Integer);
        procedure SetChristDate(y, m, d: Integer);
        procedure SetPastDate(y, m, d: Integer);

		property Book : CalcBook read fBook;
		property cycle : TPlanetCycle read fCycle;
	end;
//@+node:sonam.20180328222002.1: ** Implementation
implementation
//@+node:sonam.20180328222639.1: *3* 類型聲明
Type
    IntTab14 = array [0..13] of Integer;
    IntTab6 = array [0..5] of Integer;
	TRemainderTab = array [1..6] of Integer;
//@+node:sonam.20180328222948.1: *3* 算餘表
const
    MoonRemainderTab : TRemainderTab = (0, 50, 50, 55, 0, 45); // 第五項《脩訂再明炬》還不知
//@+node:sonam.20180328223102.1: *3* 基本方法
constructor TSysElement.Create(bk: CalcBook);
begin
	fBook := bk;
    fCycle := TPlanetCycle.Create();
    LeapMonth := false;  
end;

procedure TSysElement.SetYear(y: Integer);
begin
	fYear := y;
end;

procedure TSysElement.SetMonth(m: Integer);
begin
	fMonth := m;
end;

procedure TSysElement.SetDay(d: Integer);
begin
	fDay := d;
end;

procedure TSysElement.SetLeapMonth(l: Boolean);
begin
	fLeapMonth := l;
end;

procedure TSysElement.SetDate(y, m, d: Integer);
begin
	fYear := y;
	fMonth := m;
	fDay := d;
end;

procedure TSysElement.SetChristDate(y, m, d: Integer);
begin
    if y = 0 then
      raise EZeroYearException.Create('Zero year');

    if y > 0 then
      fYear := y + ChristianBeginYear
	else
      fYear := y + ChristianBeginYear + 1;

	fMonth := m;
	fDay := d;
end;

procedure TSysElement.SetPastDate(y, m, d: Integer);
begin
	fYear := y + ManiBeginYear;
	fMonth := m;
	fDay := d;
end;

function TSysElement.PlanetAccumMonth(y : Integer) : Int64;
begin
    result := PlanetCycle.PlanetAccumMonth(y, month, LeapMonth);
end;
//@+node:sonam.20180318215611.1: *3* 積月
//@+at 
// 所謂積月即從歷元開始加上閏月后已經經過的月數。《總滙如意寶》中採用的歷元是丁卯西元1987年，這樣計算出的積月就与歷元相關，並且僅對於歷元之后的時間才有意義。這裏採用三日共週年來計算積月，可以保證計算結果對於相當長的時期里都有意義。
// 
// 此處積月的計算采用《總滙如意寶》中第58頁中的計算方法：
// 
// <積月> ＝ [ (<曜三日共週年> * 12 + <已過月數>) * 67 + <算餘> ] div 65
// 
// 其中的<已過月數>對于藏歷三月即是零，四月是一，以此類推。
// 
//@@c
function TSysElement.AccumMonth(): Integer;
Var
    y : Integer;
begin
	y := cycle.YearInCycle(book, Week, year);  // 曜三种日共週年
    result := PlanetAccumMonth(y);
end;
//@+node:sonam.20180318215631.1: *3* 閏餘
//@+at 
// 前面計算積月時除以65之后的餘數叫作閏餘，它是用來檢驗閏月的標準：
//     当閏餘爲48或者49時，該月爲閏月
//     在閏餘爲0或者1時，積月應當減去1
//@@c
function TSysElement.MonthRemainder(): Integer;
Var
    y : Integer;
    p : Integer;
begin
	y := cycle.YearInCycle(book, Week, year);  // 曜三种日共週年
    p := (y * 12 + month) * 67;
    if (p mod 65) <= 1 then
        if not LeapMonth then
            result := (p mod 65) + 65
        else
            result := p mod 65
    else
        result := p mod 65;
end;

function TSysElement.IsExcessMonth() : Boolean;
Var
    monRemainder : Integer;
begin
    monRemainder := MonthRemainder();
    result := (monRemainder = 50) or (monRemainder = 51);
end;
//@+node:sonam.20180318215656.1: *3* 曜基數
//@+at 所謂曜基數即所求月的平朔時刻和平朔日的週日。亦即所求月之前一月三十日的值日曜（星期序數）和該太隂日結束時刻，也即平朔時刻。
// 
// 根據《歷算总汇如意宝》58頁至59頁，曜基數的計算方法如下：
// 
// <曜積月> = ((<曜三种日共週年> + <自丁卯(1987)已過週年>) * 12 + <已過月數>)
//              * 67 div 65;
// 
// <被除數> = <曜積月> mod <太隂月週期39592> * <週期數8657> mod <太隂月週期>
//  * <第一位分母7>
// 
// <曜基數> = TWeekData(被除數,0,0,0,0,0).Divide(<太隂月週期>)
// 
//@@c
function TSysElement.PlanetBase() : TWeekData;
Var
    dividend : Integer;
begin
	result := TWeekData.Create([0,0,0,0,0,0]);
	dividend := cycle.GetDividend(Week, AccumMonth()) * result.wheel[1];
	result.data[1] := dividend;
	result.Divide(cycle.MonthCycle(Week));
end;
//@+node:sonam.20180318215712.1: *3* 整零
//@+at 整數表示月亮在所求月平朔距遠地點的整日數
// 
// 根據《總滙如意寶》59頁至60頁的計算方法
// 
// <被除數> = <整零之太隂積月> mod <整零週期> * <整零週期數> mod <整零週期> * 28
// 
// <整零> = TMoonData(<被除數>,0).Divide(整零太隂月週期)
// 
// 其中，
// 
// <整零之太隂積月> = [(<自丁卯已過週年>+<整零之三日共週年>) * 12 + 已過月數] * 67 div 65
// 
//@@c
function TSysElement.IntFraction() : TMoonData;
var
    yearInCycle : Integer;
    dividend : Integer;
begin
	result := TMoonData.Create([0, 0]);
	yearInCycle := cycle.YearInCycle(book, Moon, year);
	dividend := cycle.GetDividend(Moon, PlanetAccumMonth(yearInCycle)) * result.wheel[1];
	result.data[1] := dividend;
	result.Divide(cycle.MonthCycle(Moon));
end;
//@+node:sonam.20180318215723.1: *3* 太陽基數
//@+at 太陽基數爲所求月平朔時太陽距春分點的弧長，即該時太陽的平黃經
// 
// 根據《總滙如意寶》60頁的計算方法
// 
// <被除數> = <太陽太隂積月> mod <太陽三日週期> * <太陽週期數>
//             mod <太陽三日週期> * <第一位分母>
// 
// <太陽基數> = TSunData(<被除數>,0,0,0,0).Divide(<太陽三日週期>)
// 
//@@c
function TSysElement.SunBase() : TSunData;
var
	yearInCycle : Integer;
    dividend : Integer;
begin
	result := TSunData.Create([0,0,0,0,0]);
	yearInCycle := cycle.YearInCycle(book, Sun, year);
	dividend := cycle.GetDividend(Sun, PlanetAccumMonth(yearInCycle)) * result.wheel[1];
	result.data[1] := dividend;
	result.Divide(cycle.MonthCycle(Sun));
end;
//@+node:sonam.20180318215745.1: *3* 中曜和中日
//@+at 中曜是所求日的曜次和太隂日結束的時刻，因尚未作月運行快慢的脩正，這個數值還不十分準確，隻是求定曜的一個中間步驟，故名中曜。中日爲所求的太隂日結束時太陽的平黃經，因尚未作月運行快慢的脩正，這個數值還不十分準確，隻是求定日的一個中間步驟，故名中日。
// 
// 根據《總滙如意寶》61頁：
// 
// <中曜> = <曜基數> + <曜運行速度0,59,3,4,1,365> * 日期;
// 
// <中日> = <太陽基數> + <太陽運行速度0,4,21,5,43> * 日期;
// 
//@@c
function TSysElement.PlanetMid() : TWeekData;
begin
	result := TWeekData.Create([0,59,3,4,1,365]);
	result.Multiply(day);
	result.Add(PlanetBase());
	result.Calibrate();
end;

function TSysElement.SunMid() : TSunData;
begin
	result := TSunData.Create([0,4,21,5,43]);
	result.Multiply(day);
	result.Add(SunBase());
	result.Calibrate();
end;
//@+node:sonam.20180318220223.1: *3* 月離步數表
//@+at 月離步數表是用來刻畫月亮運行遲速變化的，一個近點月大約是二十八天，月離表刻畫其中的半個週期，另外半個週期與此對稱。月離表的第一列是檢步序數，第二列叫作損益率，即諸曜從遠地點(或近日點)開始每行一步，相對於平均速度而言，超過或不及的度數，單位是弧刻。月亮的運行速度在五十四和六十四弧刻之間變化，所以與平均行度之差最大爲五。第三列叫作盈縮積，是從遠地點開始到某曜所在的宮(宿)之間，各步的損益率纍積之和。
// 
//@@c
function GetMoonMultiplier(step : Integer) : Integer;
Const
    MoonMultiplierTab : IntTab14 =
        ( -5, 5, 5, 5, 4, 3, 2, 1, -1, -2, -3, -4, -5, -5 );
Var
    i : Integer;
begin
    if step = 13 then
        i := 0
    else
        i := step+1;
    result := MoonMultiplierTab[i];
end;

function GetMoonAccum(step : Integer) : Integer;
Const
    MoonAccumTab : IntTab14 = 
        ( 0, 5, 10, 15, 19, 22, 24, 25, 24, 22, 19, 15, 10, 5 );
begin
    result := MoonAccumTab[step];
end;
//@+node:sonam.20180318220234.1: *3* 月步
//@+at 根據《總滙如意寶》中的計算方法，
// 
// <月步> = <盈縮積> + <淨行度>
// 
// 其中：
// <淨行度> = TWeekData(0,<被除數>,0,0,0,0).Divide
//                     (TMoonData.wheel[2]也即126)
// 
// 被除數的計算方法是：
// 
// <商餘> = (<整數> + <日期>) div 14
// 
// 以<商餘>查月步表之<盈縮積>，以其爲第二位建立曜數據之六位分數，以查月步表之<損益率>乘以<零數>后的得到<被除數>
// 
//@@c
function TSysElement.MoonStep() : TWeekData;
Var
	intfrac : TMoonData;
	quotient : Integer;
	remainder : Integer;
	divider : Integer;
	multiplier : Integer;
	planetMove : TWeekData;
	backStep : Boolean;
begin
	intfrac := IntFraction();
	quotient := (intfrac.data[1] + day) div 14;
	MoonBalance := not Odd(quotient);
	remainder := (intfrac.data[1] + day) mod 14;
	multiplier := GetMoonMultiplier(remainder);
	divider := multiplier * intfrac.data[2];
	if divider < 0 then
	begin
		divider := -divider;
		backStep := true;
	end
	else
		backStep := false;
	planetMove := TWeekData.Create([0,divider,0,0,0,0]);
	planetMove.Divide(126);
	result := TWeekData.Create([0,GetMoonAccum(remainder),0,0,0,0]);
	if backStep then
		result.Sub(planetMove)
	else
		result.Add(planetMove);
end;
//@+node:sonam.20180318220655.1: *3* 半定曜
//@+at 所謂半定曜，即中曜進行了月亮運行快慢脩正後的結果，由於尚未進行太陽運行快慢的脩正，還不是最后的定曜，故名半定曜
// 
//@@c
function TSysElement.SemiAccuratePlanet() : TWeekData;
Var
    mnStep : TWeekData;
begin
    result := PlanetMid();
    mnStep := MoonStep();
    if MoonBalance then
        result.Add(mnStep)
    else
        result.Sub(mnStep);
    result.Calibrate();
end;
//@+node:sonam.20180318220710.1: *3* 日躔步數表
//@+at 日躔步數表是用來刻畫太陽運行遲速變化的，太陽運行一週要經過黃道十二宮，日躔表刻畫其中的半個週期，另外半個週期與此對稱。日躔表的第一列是引檢步數，第二列叫作損益率，第三列叫作盈縮積，是從遠地點開始到某曜所在的宮(宿)之間，各步的損益率纍積之和。
// 
//@@c
function GetSunMultiplier(step : Integer) : Integer;
Const
    SunMultiplierTab : IntTab6 =
        ( -6, 6, 4, 1, -1, -4 );
Var
    i : Integer;
begin
    if step = 6-1 then
        i := 0
    else
        i := step+1;
    result := SunMultiplierTab[i];
end;

function GetSunAccum(step : Integer) : Integer;
Const
    SunAccumTab : IntTab6 = 
        ( 0, 6, 10, 11, 10, 6 );
begin
    result := SunAccumTab[step];
end;
//@+node:sonam.20180318220736.1: *3* 日步
//@+at 日步，即太陽盈縮定步，是太陽實際上比平均行度多行或者少行的弧長，
// 根據《總滙如意寶》第60頁：
// 
// <日步> = <盈縮積> + <太陽淨行>
// 
// <太陽距遠地點弧長> = <中日> - (6,45,0,0,0)
// 如果結果大於半週,則要再減去(13,30,0,0,0)
// 
// <太陽距遠地點弧長>化成弧刻，除以135化成宮數，以得到的商作為宮數查
// 日躔步數表，得到作為太陽數據第二位的<盈縮積>，以除以135的餘數作為
// 太陽數據之第二位，後面的位取<太陽距遠地點弧長>的分以下各位，乘以查表
// 得到的損益率，進位，然后每一位除以135，得到<太陽淨行>
// 
//@@c
function TSysElement.SunStep() : TSunData;
Var
	sunMiddle : TSunData;
	halfCycle : TSunData;
	step : Integer;
	remainder : Integer;
	backStep : Boolean;
	multiplier : Integer;
	divider : TSunData;
begin
	sunMiddle := SunMid();
	sunMiddle.Sub(TSunData.Create([6,45,0,0,0]));
	halfCycle := TSunData.Create([13,30,0,0,0]);
	if sunMiddle.GreatEqual(halfCycle) then
	begin
		SunSecondHalf := true;
		sunMiddle.Sub(halfCycle)
	end
	else
		SunSecondHalf := false;

	step := (SunMiddle.data[1] * 60 + SunMiddle.data[2]) div 135;
	remainder := (SunMiddle.data[1] * 60 + SunMiddle.data[2]) mod 135;

	backStep := false;
	multiplier := GetSunMultiplier(step);
	if multiplier < 0 then
	begin
		multiplier := -multiplier;
		backStep := true;
	end;
	divider := sunMiddle;
	divider.data[1] := 0;
	divider.data[2] := remainder;
	divider.Multiply(multiplier);
	divider.Divide(135);
    divider.Calibrate();

	result := TSunData.Create([0,GetSunAccum(step),0,0,0]);
	if not backStep then
        result.Add(divider)
    else
        result.Sub(divider);
end;
//@+node:sonam.20180318220750.1: *3* 定曜和定日
//@+at 根據《總滙如意寶》第68頁：
// 
// <定曜> = <半定曜> +/- <日步>
// 
// <定日> = <中日> +/- <日步>
// 
// 其中的加減是看太陽運行到了前半個週期還是後半個週期，後半週期加，前半週期減
// 
//@@c
function TSysElement.AccuratePlanet() : TWeekData; // 定曜
Var
    planetSunStep : TWeekData;
    curSunStep : TSunData;
    i : Integer;
begin
    result := SemiAccuratePlanet();
    curSunStep := SunStep();
    planetSunStep := TWeekData.Create([0,0,0,0,0,0]);
    for i := 1 to 5 do
        planetSunStep.data[i] := curSunStep.data[i];

    if SunSecondHalf then
        result.Add(planetSunStep)
    else
        result.Sub(planetSunStep);
    result.Calibrate();
end;

function TSysElement.AccurateSun() : TSunData;    // 定日
begin
    result := SunMid();
    if SunSecondHalf then
        result.Add(SunStep())
    else
        result.Sub(SunStep());
    result.Calibrate();
end;
//@+node:sonam.20180318220802.1: *3* 太陽日月宿
//@+at 太隂日月宿即太隂日結束時月亮所在之宿，亦即此時月亮的黃經
//@@c
function TSysElement.MoonStar() : TSunData;
Var
    moonDelta : Integer;
begin
    moonDelta := 54 * day;
    result := TSunData.Create([moonDelta div 60, moonDelta mod 60, 0,0,0]);
    result.Add(AccurateSun());
    result.Calibrate();
end;

// 太陽日月宿又稱伴曜月宿，是太陽日開始時月亮所在之宿和該宿中已運行的弧長，
// 因爲值日曜的起訖是按太陽日計算的，與太陽日共始終的，所以
// “伴曜”就是按太陽日計算的意思，亦即天明時月亮的黃經。
function TSysElement.SonarMoonStar() : TWeekData;
Var
    curMoon : TSunData;
    planet : TWeekData;
    i : Integer;
begin
    curMoon := MoonStar();
    result := TWeekData.Create([0,0,0,0,0,0]);
    for i:=1 to 5 do
        result.data[i] := curMoon.data[i];
    planet := AccuratePlanet();
    planet.data[1] := 0;
    result.Sub(planet);
end;
//@+node:sonam.20180318220831.1: *3* 會合
//@+at 根據《總滙如意寶》
// 
// <會合> = <太陽日月宿> + <定日>
// 
//@@c
function TSysElement.Intersect() : TSunData;
Var
    curSonarMoonStar : TWeekData;
    i : Integer;
begin
    result := TSunData.Create([0,0,0,0,0]);
    curSonarMoonStar := SonarMoonStar();
    for i:=1 to 5 do
        result.data[i] := curSonarMoonStar.data[i];

    result.Add(AccurateSun());
    result.Calibrate();
end;
//@+node:sonam.20180318220840.1: *3* 作用
//@+at 
// 作用及其前后分是把所推算的這個太隂日分成兩半，每一半是一個作用的數值。其算法是：用前一個太隂日減去六十，差數是在前一個太隂日里已走過的，這個差數與今天的太隂日數值相加，除以二，商數即這個時間的作用的數值。前分是前一個太隂日終了後還要經過這么長的時間，后分表示這個數值終了後的作用的數值。
// 
//@@c
function TSysElement.Effect() : Integer;
begin
    result := (day * 2 - 1) mod 7;
end;
//@-others

end.
//@-leo
